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INTRODUCTION 


Analysis  and  detailed  understanding  of  the  flowfield  produced  by  a 
propeller/body/appendage  combination  is  of  interest  in  a  number  of  practical 
applications.  Examples  Include  the  influence  of  the  propeller  on  the  body 
pressure  distribution  to  estimate  the  thrust  deduction,  the  effects  of  non- 
uniform  inflow  on  propeller  performance,  the  prediction  of  the  near-wake 
profiles,  cyclic  loading  to  produce  vibrations  and  unsteady  bearing  forces 
and  the  influence  of  a  rudder  or  other  appendages  on  the  propeller  and  vice 
versa.  We  have  been  conducting  a  research  plan  that  addresses  two  important 
parts  of  the  general  problem:  1)  numerical  analysis  of  the  3~D,  turbulent 
flow  produced  by  propeller/ body/appendage  combinations  and  2)  studies  of  the 
large-scale  .and  unsteady  structure  of  3~D,  turbulent  boundary  layer  flows 
produced  by  hull/appendage  junctions.  This  section  contains  some  background 
material  on  these  research  areas. 

Numerical  Solution  of  3~Dt  Turbulent  Propeller  Flows 

We  consider  first  the  current  state  of  knowledge  for  analysis  of  3~D, 
turbulent  flows  near  and  through  propellers.  Until  quite  recently,  only 
approximate  treatments,  e.g..  Refs.  (1)  -  (5)  were  available.  These 
analyses  involved  one  or  more  of  the  following  restrictive  assumptions:  the 
flow  was  assumed  inviscid;  the  propeller  was  represented  as  an  actuator  disk 
with  constant  thrust  over  the  disk;  the  flow  was  taken  as  laminar  and/or  the 
effects  of  the  propeller  on  the  flowfield  were  assumed  small  enough  to 
permit  a  linearization  of  the  equations  of  motion. 

The  numerical  solution  procedure  reported  in  Refs.  (6)  and  (7)  had  as 
its  goal  the  development  of  a  realistic  treatment,  holding  simplifying 
assumptions  and  approximations  to  a  minimum.  The  work  was  based  on  the 


full,  axisymmetric,  mean  (in  the  turbulence  sense),  unsteady  Navier-Stokes 
equations.  In  order  to  place  some  bounds  on  the  scope  of  the  effort, 
however,  some  simplifications  were  necessary.  The  first  was  the  assumption 
of  an  actuator  disk  model  for  the  propeller,  although  arbitrary  radial 
variations  of  thrust  were  allowed.  Second,  the  flow  was  taken  as 
axisymmetric.  Third,  turbulent  transport  processes  were  described  by  an 
integrated,  turbulence-kinetic-energy  (TKE)  equation,  which  was  used  to 
predict  an  eddy  viscosity  distribution.  The  eddy  viscosity  was  allowed  to 
vary  only  in  the  streamwise  direction.  The  simulation  of  a  given  propeller 
for  the  purposes  of  a  viscous  numerical  calculation  was  considered  in  detail 
for  the  first  time  in  this  work.  The  unsteady  equations  of  motion  were  cast 
in  terms  of  a  stream  function,  one  vorticity  component,  and  the  peripheral 
velocity.  The  vorticity  equation  and  the  peripheral  momentum  equation  were 
solved  by  an  Alternating  Difference  Implicit  technique,  and  the  Poisson 
equation  for  the  stream  function  was  solved  by  Direct  Matrix  Reduction. 

Comparison  of  the  predictions  of  that  procedure  with  laboratory  data 
for  an  axisymmetric,  slender,  upstream  body  with  no  appendages  showed  good 
agreement  for  the  axial  velocity  profile  at  X/D  -  2.  The  swirl  velocity 
component  was  underpredicted. 

From  the  point  of  view  of  direct  applications  to  problems  of  practical 
interest,  the  biggest  limitation  to  the  analysis  described  above  is  that  it 
was  restricted  to  two-dimensional,  axisymmetric  flows.  Actual  propeller- 
driven  vessels  have  either  a  three-dimensional  body  ahead  of  the  propeller 
and/or  appendages  such  as  rudders  or  control  planes  near  the  propeller  that 
render  the  flow  basically  three-dimensional,  even  if  the  propeller  is  still 
taken  as  an  actuator  disk,  neglecting  the  explicit  influence  of  individual 
blades.  For  many  cases,  the  assumption  of  an  actuator  disk  representation 
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of  the  propeller  remains  appropriate  while  the  restriction  to  axisymmetric 
flows  does  not. 

The  research  program  undertaken  in  1983  under  ONR  support  is  a  step-by- 
step  approach  to  the  development  of  a  method  for  handling  propeller 
flowfields  with  three-dimensional  inflows  and  the  effects  of  appendages  in 
situations  of  increasing  complexity  and  practical  realism,  while  retaining  a 
general  actuator  disk  model  of  the  propeller.  The  problem  treated  in  the 
first  year  was  an  isolated  propeller  operating  in  an  approach  flow  with  a 
linear  velocity  gradient.  The  numerical  solution  procedure  was  based  on  the 
Finite  Element  Method  (FEM)  for  the  steady  Navier-Stokes  equations  written 
in  primitive  variables  following  the  suggestions  of  Engelman  (Refs.  (8),  (9) 
and  (10)).  Turbulence  modeling  was  through  an  integrated  form  of  the 
Turbulent  Kinetic  Energy  (TKE)  equation.  In  order  to  have  flowfield 
measurements  against  which  the  adequacy  of  the  numerical  predictions  could 
be  tested,  a  series  of  wind  tunnel  experiments  was  run. 

Calculations  have  been  run  for  two  cases  where  experimental  data  are 
available  -  one  with  a  uniform  approach  velocity  as  a  baseline  case  and  one 
with  the  linear  velocity  gradient  approach  velocity  profile.  Comparisons  of 
predictions  and  measurements  are  given  in  Ref.  (11)  (see  the  report  in  the 
Appendix)  for  streamwlse  velocity  profiles  just  behind  the  propeller  disk 
(X/D  =  .025).  The  numerical  (FEM)  predictions  give  a  good  description  of 
the  flow. 

3-D  Turbulent  Flows  from  Hull  Appendage  Junctions 

When  a  laminar  or  turbulent  boundary  layer  on  a  surface  encounters  a 
protuberance  of  strut  projecting  from  that  surface,  a  "horseshoe"  vortex  is 
formed  at  the  junction  of  the  two  surfaces.  Upstream  of  the  protuberance 


4 


the  time-averaged  vorticity  direction  in  the  boundary  layer  is  spanwise 
across  the  surface.  In  order  to  satisfy  the  vortex  theorems  of  fluid 
dynamics,  streamwise  legs  of  this  upstream  vorticity  stretch  around  the 
protuberance  in  a  horseshoe  shape  with  each  leg  having  vorticity  of  opposite 
rotational  sense.  This  type  of  three-dimensional  turbulent  boundary  layer 
occurs  in  a  number  of  practical  cases,  such  as  in  turbomachinery  blade  and 
end  wall  flows,  aircraft  wing  and  body  junction  flows,  and  ship  appendage 
and  hull  Junction  flow. 

Dean  (Ref.  12)  and  Humphreys  and  vandenBerg  (Ref.  13)  reviewed 
available  data  for  this  type  of  flow  in  order  to  select  the  best  test  cases 
for  evaluating  computational  methods  (Ref.  14).  Dean  pointed  out  that 
little  turbulence  data  had  been  obtained  during  50  years  of  measurements. 
Consequently,  only  the  Shabaka  (Ref.  15)  flow  around  a  flat  plate  with  a 
semi-elliptical  nose  was  recommended  as  a  computational  test  case  for  the 
secondary  flow  produced  by  the  streamwise  vortex  legs  well  downstream  of  the 
nose.  Humphreys  and  vandenBerg  recommended  the  Dechow  (Ref.  16)  flow  as  a 
test  case,  in  which  turbulence  measurements  were  made  in  front  of  and  beside 
the  protuberance  or  obstruction,  but  not  in  the  presence  of  adverse  pressure 
gradients.  In  this  case  the  protuberance  was  a  320  mm  diameter  semi¬ 
circular  cylinder  with  a  streamlined  afterbody  that  prevented  premature 
separation  upstream  of  the  trailing  edge.  Pierce’s  group  (Ref.  17)  have 
made  some  time-averaged  measurements  around  the  nose  of  a  circular-nosed 
strut. 

From  a  computational  viewpoint,  the  nose  region  is  the  easiest  to 
calculate  since  the  Incident  boundary  layer  mean  vorticity  distribution  and 
an  inviscid  rotational  solution  technique  largely  determine  the  rapidly 
distorted  flowfield  around  the  nose.  Computational  methods  for  the  mean 


flow  in  the  corners  around  a  protuberance  are  not  in  a  satisfactory  state, 
even  for  the  zero  pressure  gradient  case  (Ref.  14).  Dean  pointed  out  the 
need  for  measurements  of  the  streamwise  vortex  corner  flow  in  the  presence 
of  adverse  pressure  gradients  for  increased  understanding  and  for 
improvement  and  verification  of  computational  methods.  In  most  practical 
applications  these  adverse  pressure  gradients  cause  the  vortex  strength  to 
be  attenuated  and  the  vortex  to  wander  away  from  the  corner  (Ref.  19). 

In  other  work  on  streamlined  struts.  Barber  (Ref.  20)  showed  that  the 
size  of  the  separated  flow  zone  was  inversely  related  to  the  incident 
boundary  layer  thickness.  Mehta  (unpublished,  but  referenced  in  Ref.  19) 
showed  that  the  streamwise  corner  vortex  diameter  and  distance  from  the  body 
were  dependent  on  the  leading-edge  shape.  Recently,  Oguz  (Ref.  21)  made 
measurements  on  the  flow  around  a  body  of  constant  thickness,  a  case  similar 
to  Shabaka's. 

In  all  of  these  studies  only  time-averaged  measurements  were  made.  No 
quantitative  measurements  were  made  of  flow  regions  that  contain  adverse 
pressure  gradients  and/or  separation.  No  spectra  of  the  turbulence  energy 
containing  frequencies  were  obtained.  No  examination  has  been  made  of  the 
organized  time-varying  coherent  structure  of  this  type  of  turbulent  boundary 
layer,  such  as  has  been  done  for  mean  two-dimensional  turbulent  boundary 
layers. 

Recently,  Rood  at  DWTNSRDC  reported  relatively  narrow  frequency 
bandwidth  turbulent  fluctuations  in  the  three-dimensional  boundary  layer 
downstream  of  a  streamlined  protuberance.  In  some  cases  this  downstream 
boundary  layer  passes  over  other  protuberances  and  structural  elements  and 
Interacts  with  them  in  an  undesirable  and  as  yet  unpredictable  way.  Noise 
is  generated  when  large-scale  turbulent  structures  in  this  boundary  layer 
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interact  with  these  structural  elements.  When  a  rotor  blade  encounters  the 
wake  of  a  strut,  tonal  noise  is  generated.  It  is  thought  that  the  organized 
time-varying  coherent  structure  of  these  flows  may  be  able  to  explain  such 
undesirable  behavior  and  to  aid  in  future  optimal  designs  of  these 
structural  elements. 

A  time-averaged  description  of  the  three-dimensional  flow  downstream  of 
a  protuberance,  whether  from  experiments  or  a  calculation  method,  provides 
no  information  on  the  large-scale  time-dependent  motions.  Since  these  large- 
scale  time-dependent  motions  would  interact  with  downstream  moving  rotor 
blades  differently  than  for  the  mean  flow,  information  is  needed  on  the 
nature  of  these  large-scale  structures  of  a  three-dimensional  flow 
downstream  of  a  protuberance. 

Thus,  although  there  have  been  previous  time-averaged  measurements  of 
this  type  of  flow  in  flow  regions  without  pressure  gradients  or  separation, 
there  has  been  no  investigation  of  the  role  of  the  time-dependent  organized 
turbulent  flow  structure,  and  no  quantitative  measurements  of  flow  regions 
that  contain  adverse  pressure  gradients  and/or  separation. 

At  this  time  we  can  only  conjecture  about  how  a  protuberance  affects 
the  organized  turbulent  flow  behavior,  based  upon  our  experience  from  mean 
two-dimensional  boundary  layers.  Large-scale  coherent  structures  from  the 
upstream  boundary  layer  will  be  stretched  around  the  protuberance. 
Depending  upon  the  incident  Reynolds  number  and  the  geometry  of  the 
protuberance,  certain  size  or  frequency  structures  may  grow  or  be  attenuated 
and  a  more  narrow  bandwidth  spectrum  may  result.  Experiments  over  a  range 
of  conditions  are  needed  to  determine  when  narrow  band  spectra  occur. 

Adverse  pressure  gradients  on  the  downstream  side  of  a  protuberance  or 


strut  may  cause  the  individual  structures  to  move  away  from  the  corner  and 
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agglomerate  by  pairing,  since  Shabaka  and  Bradshaw  (Ref.  19)  have  noted  that 
the  mean  flow  vortex  does  wander  from  the  corner  under  these  conditions. 
Separation  is  likely,  with  rapidly  growing  large-scale  structures  supplying 
reversed  flow,  which  also  occurs  in  the  mean  two-dimensional  separation 
case.  In  fact  there  are  rw  detailed  mean  flow  turbulence  measurements 
available  for  this  separated  flow  in  the  corner. 

The  flow  downstream  of  the  protuberance  or  strut  is  a  three-dimensional 
turbulent  boundary  layer.  If  the  downstream  flow  is  subjected  to  an  adverse 
pressure  gradient  as  is  often  the  case,  the  protuberance  corner  separation 
may  not  reattach,  but  instead  the  separation  can  spread  laterally  across  the 
boundary  layer.  Even  if  the  separated  corner  flow  does  reattach,  this 
momentum  and  energy  deficient  region  of  the  flow  is  more  prone  to  detachment 
when  subjected  to  an  adverse  pressure  gradient  downstream.  In  cases  where 
these  types  of  flow  interact  with  another  structural  element  downstream,  the 
flow  structure  is  clearly  important  in  designing  the  structural  element  for 
optimal  performance. 

In  the  current  research  program,  which  began  May  15,  1984,  measurements 
are  and  will  be  made  on  turbulent  boundary  layers  around  protuberances  or 
struts  and  downstream,  including  strong  adverse  pressure  gradient  and 
detached  flow  regions.  Information  on  the  temporal  flow  behavior  will  be 
obtained  by  means  of  spectral  analysis  and  space-time  correlations  of 
surface  pressure  transducer  and  hot-wire  anemometry.  Time-averaged 

measurements  of  U,  V,  W,  -uv,  -uw,  u* ,  v^  and  w*  will  be  obtained  in  detail 
for  the  region  around  the  strut  and  the  downstream  three-dimensional 
turbulent  boundary  layer,  including  zones  of  detached  flow.  Surface 
pressure  spectra  can  be  used  to  estimate  sound  pressure  levels  outside  of 


the  turbulent  flow. 


VELOCITY  COMPONENT,  VY/VT,  AT  Y=0.  RPM=UOO  ,X/DP=.025 


Figure  7. 


Swirl  Velocity  Profile  Behind  the 
Propeller  at  x/d=.C25  for  a  Vertical  Trav¬ 
erse  . 
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VELOCITY  COMPONENT,  VY/VF,  AT  Z=0,  RPM=1400  ,X/DP=.025 


Figure  5.  Radial  Velocity  Profile  Behind  the 
Propeller  at  x/d=.025  for  a  Horizontal 
Traverse . 
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VELOCITY  COMPONENT,  VX/VF,  AT  Z=0,  RPM=1400  ,X/DP=.025 
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wire  probes  also  show  normal  zero-pressure-gradient  flow  behavior  without 
any  preferred  frequency  (Ref.  2*1). 

In  work  sponsored  by  NASA,  McGrath  and  Simpson  (Ref.  23)  have  obtained 
surface  pressure  spectra  (Figures  14  and  15),  pressure  fluctuation  wave 
speeds,  and  streamwise  and  spanwise  coherence  measurements.  Figure  14  shows 
that  the  surface  pressure  spectra  at  lower  frequencies  agree  with  one 
another  in  terms  of  the  outer  region  variables  non-dimenslonalizing 
variables  and  the  displacement  thickness  6.  Figure  15  shows  that  the 

higher  frequency  surface  pressure  fluctuation  spectra  correlate  when  non- 
dimensionalized  on  wall  variables  U  and  v/U  .  No  discrete  preferred 

T  T 

frequency  has  been  observed  in  any  of  these  surface  pressure  fluctuation 
data. 

Figure  16  shows  the  shape  of  the  cylindrical  appendage  that  is 
currently  being  fabricated  on  a  computer-controlled  milling  machine  and  will 
be  used  in  the  wind  tunnel  shown  in  Figure  17.  The  nose  of  this  12  inch 
chord  protuberance  is  a  1.5  to  1.0  half  ellipse  while  the  tall  portion 
beginning  at  the  maximum  thickness  has  the  shape  of  a  NACA  0020  airfoil. 
The  2.83  inches  maximum  thickness  protuberance  produces  some  blockage 
effects.  In  particular,  the  flow  in  the  trailing  portions  of  the  corners  of 
this  Junction  will  not  be  subjected  to  the  same  streamwise  pressure 
gradients  as  in  an  infinite  stream  unless  blockage  effects  are  properly 
accounted  for.  The  wall  Jet  and  suction  boundary  layer  control  system  will 
be  used  to  minimize  blockage  effects.  As  shown  in  Figure  18,  inviscid  flow 
streamlines  for  this  body  in  an  infinite  stream  show  that  there  is  some 
streamline  curvature  18  inches  from  the  body.  The  boundary  layer  wall  Jet 
and  suction  system  will  be  used  to  suck  off  flow  and  produce  this  streamline 
curvature.  Mean  surface  pressures  calculated  from  potential  flow  will  be 


than  30  boundary  layer  thicknesses  downstream  of  location  A  to  insure  that  a 


normal  two-dimensional  zero-pressure-gradient  boundary  layer  exists  just 
upstream  of  the  protuberance. 

Downstream  of  the  protuberance  the  flow  will  reattach  to  some  degree. 
Downstream  of  location  B  the  flow  is  subjected  to  a  strong  adverse  pressure 
gradient  that  can  be  adjusted  by  the  upper  wall.  The  momentum  deficient 
wake  from  the  protuberance  will  be  the  first  part  of  the  flow  to  detach  in  a 
mean  three-dimensional  pattern.  Judging  from  the  results  from  mean  two- 
dimensional  flow  the  detached  flow  structure  is  dominated  by  the  large-scale 
structures. 

+ 

Figure  1  1  shows  that  upstream  of  the  protuberance  the  classical  U  vs. 

+ 

Y  "law-of-the-wall"  velocity  profile  exists  near  the  wall.  The  surface 
skin  friction  factor  obeys  the  classical  Ludwieg-Tillmann  correlation  (Ref. 
22).  The  Reynolds  shearing  stress  distributions  across  this  boundary  layer 
were  obtained  by  cross-wire  and  triple-wire  hot-wire  anemometer  probes  (Ref. 
24)  and  are  in  good  agreement  with  earlier  distributions  for  this  type  of 

flow  (Figure  12).  Measurements  of  u*  by  single,  cross-wire,  and  triple-wire 
hot-wire  anemometers  are  in  good  agreement  as  discussed  by  Ahn  (Ref.  22). 

Measurements  of  v*  by  cross-wire  and  triple-wire  probes  agree  away  from  the 

wall  where  triple-wire  measurements  are  valid  (Ref.  24).  Measurements  of  w* 
are  in  agreement  with  earlier  turbulent  boundary  layer  results. 

Figure  13  shows  that  u  fluctuation  spectra  from  single  hot-wire  data 
have  normal  wideband  features  without  any  preferred  single  frequency  (Ref. 
22).  This  Indicates  no  appreciable  single  frequency  unsteadiness  in  this 


wind  tunnel  flow.  Spectra  for  v,  -uv  and  w  obtained  by  cross  and  triple 
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3-D  Turbulent  Flows  from  Hull  Appendage  Jcinctiona 

During  the  first  phase  of  this  work,  our  goal  is  to  determine  the  time- 
dependent  features  of  protuberance  flows  for  a  range  of  practical 
conditions.  Turbulence  spectra,  cross-spectra,  and  space-time  correlations 
of  hot-wire  anemometer  signals  of  u  around  the  protuberance  and  downstream 
will  be  used*  to  obtain  this  information.  Hot-wire  anemometer  measurements 
upstream  of  the  protuberance  and  downstream  of  reattachment  will  be  made  to 
document  each  condition. 

During  the  past  first  year  of  this  phase  (May  15,  1984  through  May  14, 
1985)f  detailed  turbulence  measurements  have  been  made  of  the  zero-pressure 
gradient  boundary  layer  that  approaches  the  protuberance  (Refs.  22,  23  and 
24).  The  temporal  structure  of  this  type  flow  has  been  examined  over  a 
range  of  momentum  thickness  Reynolds  numbers  Re^  between  2500  and  7000. 

Spanwise  and  streamwlse  space-time  correlations  and  cross-spectral 

distributions  of  fluctuations,  u* ,  v*,  w*  and  -uv  Reynolds  stresses,  and 
surface  pressure  spectra  have  been  obtained  to  thoroughly  document  the  time- 
dependent  behavior  of  zero-pressure  gradient  turbulent  boundary  layers 
approaching  a  protuberance.  This  is  extremely  important  so  that  we  can  know 
if  the  time-dependent  narrow  frequency  bandwidth  motions  observed  by  Rood 
are  due  to  the  upstream  zero-pressure  gradient  boundary  layer  of  due  to  the 
presence  of  the  protuberance. 

Figure  9  shows  a  sideview  schematic  of  the  26  feet  long,  3  feet  wide 
wind  tunnel  test  section  used  in  these  experiments.  The  flow  enters  from 
the  left,  forming  a  mean  two-dimensional  turbulent  boundary  layer  on  the 
flat  floor  of  the  test  section.  The  upper  wall  position  is  adjusted  to 
produce  a  zero  pressure  gradient  boundary  layer  between  locations  A  and  B 
(Figure  10).  A  cylindrical  protuberance  of  interest  will  be  placed  more 


problem.  We  view  the  results  as  "preliminary",  because  they  have  only  very 
recently  been  made  and  are  still  under  study.  We  used  a  finite  element  grid 
with  12  perip'  ®ral  planes  which  means  that  the  pie-shaped  mesh  piece  in 
Figure  1  occupied  one  whole  element.  It  is  likely  that  one  or  two  extra 
peripheral  planes  should  be  added  in  that  area.  Work  to  that  end  is  now 
underway.  Despite  these  reservations,  the  results  obtained  to  date  to  be 
shown  below  are  really  quite  good,  and  we  look  for  our  final  calculations  to 
provide  a  truly  excellent  simulation  of  this  complex  flow. 

In  Figure  3.  we  show  measurements  and  calculations  for  the  axial 
velocity  profile  right  behind  the  propeller  (X/D  -  0.025)  along  a  radial 
line  in  the  horizontal  direction  in  Figure  1.  For  that  cut,  the  flow  will 
not  have  passed  through  the  extra  pie-shaped  mesh  segment  before  passing 
through  the  propeller.  Figure  4  gives  the  swirling  velocity  comporent 
profile  along  the  same  radial  cut,  and  Figure  5  shows  the  corresponding 
radial  velocity.  One  can  note  the  generally  good  agreement  except  near  the 
hub  and  drive  shaft  (Y/D  <  0.2). 

Measurements  and  predictions  along  a  vertical  cut  (see  Figure  1),  i.e. 
for  flow  that  has  passed  through  the  extra  wedge  of  mesh  before  entering  the 
propeller,  at  X/D  0.025  are  shown  in  Figures  6,  7  and  8.  The  agreement  is 
good  but  not  so  good  as  for  the  other  radial  cut.  Again,  we  attribute  this 
to  a  need  for  better  element  grid  resolution  in  this  sector. 
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Numerical  Solution  of  3~D.  Turbulent  Propeller  Flows 

In  the  original  (1983)  plan,  work  for  this  year  would  have  involved  an 
axisymmetric ,  sharp-tailed,  slender  body  ahead  of  the  propeller.  The 
approach  flow  was  to  be  taken  as  three  dimensional,  i.e.  varying  with  angle 
around  the  body,  0.  This  was  called  Problem  #2.  This  represented  a  very 
large  Jump  in  complexity  beyond  Problem  #1  -  a  propeller  in  a  linear  shear 
flow.  Our  experience  with  Problem  #1  and  consultations  with  DTNSRDC 
personnel  indicated  that  the  original  Problem  #2  should  properly  be  split 
into  two  phases  -  Problems  #2A  and  #2B.  Problem  #2A  involves  a  peripherally 
varying  inflow  velocity  profile  similar  to  that  expected  in  the  original 
Problem  #2  except  that  the  solid  body  is  removed.  The  desired  inflow 
velocity  distribution  is  generated  by  a  wire  mesh  disk  with  radial  and 
peripheral  variations  in  mesh  spacing  (see  Figure  1).  This  method  is 
commonly  used  in  propeller  testing.  Problem  #2A  was  treated  this  year. 
Calculated  results  were  compared  to  wind  tunnel  measurements  using  the  same 
propeller  as  for  Problem  #1 . 

A  photograph  of  the  test  apparatus  in  the  wind  tunnel  is  shown  in 
Figure  2.  We  are  using  the  same  propeller  as  for  our  earlier  work,  but  the 
test  conditions  are  different.  The  freestream  velocity  is  44  ft/sec  and  the 
propeller  is  nominally  "self-propelled"  cancelling  the  drag  of  the  screen 
disk.  The  rpm  was  1400  giving  J  <»  .96.  The  propeller  drive  system  and  the 
thrust  and  torque  measurements  are  improvements  over  our  previous  devices. 
We  have  measured  both  mean  flow  and  turbulence  quantities. 

To  this  point,  we  have  made  complete,  but  preliminary,  calculations 
using  the  basic  methods  described  in  Ref.  (11)  extended  to  this  complex 


Freestream  Velocity  Distributions  Along  Wind  Tunnel  Test 
Section  Centerline.  Case  A  Used  in  Protuberance-Hull 
Junction  Flow  (Ref.  22). 


Zero-Pressure-Gradient  Turbulent  Boundary 
Upstream  of  the  Appendage  (Ref.  22). 


6  14m 


Figure  12.  Reynolds  Shearing  Stress  -uv/U^^  vs.  Y/6  Profiles  from 

Cross-Wire  and  Triple-Wire  Anemometers  along  the  Test 
Section  Floor  for  Case  A.  Increasing  Streamwise  Locations 
from  Top  to  Bottom  of  Figure.  Note  Displaced  Ordinates. 
Solid  Symbols  for  Triple-Wire  Results  at  4.255  m.  (Ref.  24) 
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y/i  •  0.172 


y/S  -  0.297 


y/i  -  0.523 


y/4  •  0.925 


Figure  13.  Streamwise  Velocity  Fluctuation  Spectra  at  X  =  3.289  m  for 
Various  Y/6  Positions  within  the  Case  A  Zero-Pressure- 
Gradient  Boundary  Layer.  Solid  Straight  Lines  Denote 
-5/3 

n  Inertia  Subrange  for  These  Spectra. 


Surface  Pressure  Spectra  in  Outer  Region  Non-Dimensional 
Variables  for  the  Case  A  Zero-Pressure-Gradient  Turbulent 
Boundary  Layer  Along  the  Test  Section  Floor.  Note  Good 
Agreement  Among  Data  for  0.1  <  106, /U.  <  1  (Ref.  23). 
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Figure  15 


cuV/Ut* 


Surface  Pressure  Spectra 
Variables  for  the  Case  A 
Boundary  Layer  Along  the 
Agreement  Among  Data  for 


in  Wall  Region  Non-Dimensional 
Zero- Pressure-Gradient  Turbulent 
Test  Section  Floor.  Note  Good 
cov/U  2  >  0.1  (Ref.  23). 


Body  Shape  of  the  12  Inches  Chord  Cylindrical  Appendage  To  Be  Used  in  the 
Protuberance-Hull  Junction  Flow  Experiments.  Nose  is  1.5:1  Ellipse  While 
Trailing  Portion  from  Maximum  Thickness  is  NACA  0020  Airfoil  Shape. 


Plan  View  of  Wind  Tunnel 

Top  View 
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Figure  17.  Plan  View  of  Wind  Tunnel  Test  Section  Showing  Side  Wall  Jet 
and  Suction  Boundary  Layer  Control  Slots  for  Reducing 
Blockage  Effects  of  Appendage. 


o 

OJ 


32 


tn 


o 

■2 


in 


in 


T3 

OJ 

O)  +j  c 

C  to 
•f-  3 

C  _J  'OT)  c 
0)  "So 

jr  T3  <t 
3  C  +J 

(O  01  3 

CO  ^ 


“O 

s- 

OJ 

o 

4->' 

s. 

cn 

•r— . 

to 

3 

zs 

a><4- 

o 

•p- 

o 

E 

u. 

OJ 

OJ 

Of 

4^ 

s. 

c 

c 

to 

3 

to 

t/S 

to 

c 

s. 

OJ 

s 

OJ 

V 

o 

4-> 

.o 

Ou 

x: 

c 

oo 

OJ 

4^ 

c 

o 

c 

ns 

>1 

o 

OJ 

“O 

to 

o 

s: 

o 

4-> 

ca 

c 

s. 

■o 

OJ 

OJ 

to 

5x 

(j 

ro 

OJ 

<o 

to 

C 

s. 

I 

o 


o% 

C  i/> 
•r-  a; 

c  ^ 

o 

ns  C 


o  Q-  ^  >  <U 

•f-  <u  c  c 

(/)  L. 

c  <0  ' 

O)  ro  "O 
E  •/> 

•I-  tn  3 

O  •!-  O 

U  CQ 

’Si  ^  - - 

<  r-  c 

ns  O 
a;  3:  tt5  r— 

x:  •  s  -f  I 

4->  E  O) 

fO  “D  II 

M 

•J^ 

Qj  cn 
c  c 
•*“  o 


$-  o 
o  i- 


t/) 
to 

a;  QJ  2 

c  •*->  o 


ns  OJ 
0)  c  ^ 


C  S-. 
O 

5 

o  c 

I—  ••-  >» 

Li.  S. 

”0  w 
-o  OJ  03 
•F«  S 

u  S-  = 
to  03  >> 

•r-  s  oo 

>  -Q 
C  3  ^ 

I— •  tio  o 


00 


5  < 

0;  c 
i-  o 
•♦->  *1- 
00  •*-> 
u 

xz  <v 
u  s- 

3  -f- 

oo  o 

O)  O) 
U  to 
3  •»- 
“D  5 
O  E 
i.  ns 
a.  O) 
s. 
o 

40  OO 


O) 

S- 

3 

cr> 


o  in 

(\J 


if  o 


in 


o 


N 


REFERENCES 


Sparenberg,  J.A.,  "On  the  Potential  Theory  of  the  Interaction  of  an 
Acutator  Disk  and  a  Body,"  Journal  of  Ship  Research,  December  1972. 

Sparenberg,  J.A.,  "On  the  Linear  Theory  of  an  Actuator  Disk  in  a 
Viscous  Fluid,"  Journal  of  Ship  Research,  March  1974. 

Kim,  K-H.  and  Kobayoshi,  S.,  "Pressure  Distribution  on  Propeller  Blade 
Surface  Using  Numerical  Lifting  Line  Theory,"  SNAME  Propellers  '84 
Meeting,  May  1984. 

Huang,  T.  et  al . ,  "Propeller/Stern/Boundary-Layer  Interaction  on 
Axisymmetrlc  Bodies,"  DTNSRDC  Rept.  76-0113,  December  1976. 

Cox,  B.D.  and  Hansen,  A.G.,  "A  method  for  Predicting  Thrust  Deduction 
using  Propeller  Lifting  Surface  Theory,"  DTNSRDC  Rept.  77-0087, 
November  1977. 

Schetz,  J.A.  and  Favin,  S.,  "Numerical  Solution  for  the  Near  Wake  of  a 
Body  with  Propeller,"  Journal  of  Hydronautics,  Vol.  11,  No.  4,  October 
1977,  pp.  136-141. 

Schetz,  J.A.  and  Favin,  S.,  "Numerical  Solution  of  a  Body-Propeller 
Combination  Flow  Including  Swirl  and  Comparisons  with  Data,"  J.  of 
Hydronautics,  Vol.  13,  No.  2,  pp.  46-51,  April  1979. 

Bercovier,  M.  and  Engelman,  M.S.,  "A  Finite  Element  for  Incompressible 
Fluid  Flows,"  J.  Comp.  Physics,  30,  181  (1979). 

Engelman,  M.S.,  Sani,  R.L.  and  Gresho,  P.M.,  "The  Implementation  of 
Normal  and  Tangential  Velocity  Boundary  Conditions  in  Finite  Element 
Codes  for  Incompressible  Fluid  Flow,"  to  appear  Int.  J.  Num.  Meth. 
Fluids,  1981 . 

Engelman,  M.S.,  Strang,  G.  and  Bathe,  K.J.,  "The  Application  of  Quasi- 
Newton  Methods  in  Fluid  Mechanics,"  Int.  J.  Num.  Meth.  Eng.,  1983. 

Pelletier,  D.  and  Schetz,  J.A.,  "A  Navier-Stokes  Calculation  of  3D 
Turbulent  Flow  Near  a  Propeller  in  a  Shear  Flow,"  AIAA  Paper  85-0365, 
January  1985. 

Dean,  R.  "Turbulent  Secondary  Flow  of  the  First  Kind,"  see  Kline  et. 
al.  (1981),  pp.  139-154,  1980. 

Humphreys,  D.A.  and  vandenBerg,  B.,  "Three-Dimensional  Turbulent 
Boundary  Layers,"  see  Kline  et  al.  (1980),  pp.  162-169,  1980. 

Kline,  S.J.,  Cantwell,  D.J.  and  Lllley,  G.M.,  ed..  The  1980-1981  AFOSR- 
HTTM-Stanford  Conferences  on  Complex  Turbulent  Flows;  Comparison  of 
Computation  and  Experiment,  Vol.  1,  2  and  3,  Stanford  University  Dept. 
Mechanical  Engineering,  1981. 


34 


15.  Shabaka,  "Turbulent  Flow  in  an  Idealized  Wing-Body  Junction," 

Ph.D.  Dissertation,  Imperial  College  of  Science  and  Technology,  London, 
1979. 

16.  Dechow,  R.,  Dissertation  University  Karlsruhe,  see  Dechow,  R.  and 
Felsch,  K.O.,  Turb.  Shear  Flows  Sym.  (1977),  pp.  9.11-9.20,  Pa.  State 
Univ.  (1977). 

17.  Mcallister,  J.E.,  Pierce,  F.J.  and  Tennant,  M.H.,  "Direct  Force  Wall 
Shear  Stress  Measurement  in  Pressure-Driven  Three-Dimensional 
Turbulent  Boundary  Layers,"  J.  Fluids  Engrg.,  104,  pp.  150-155,  1982. 

18.  Rubel,  A.,  "Computational  Technique  for  the  Calculation  of  Inviscid 
Rotational  Jet  Impingement,"  AIAA  Paper  78-1212,  1978. 

19.  Shabaka,  I.M.M.A.  and  Bradshaw,  P.,  "Turbulent  Flow  Measurements  in  an 
Idealized  Wing/Body  Junction,"  AIAA  Journal,  19,  pp.  131-132,  1981. 

20.  Barber,  T.J.,  "An  Investigation  of  Strut-Wall  Intersection  Losses,"  J. 
Aircraft,  15,  pp.  576-681,  1978. 

21.  Oguz,  E.A.,  "An  Experimental  Investigation  of  the  Turbulent  Flow  in  the 
Junction  of  a  Flat  Plate  and  a  Body  of  Constant  Thickness,"  Ph.D. 
Dissertation,  Georgia  Inst.  Tech.,  1981. 

22.  Ahn,  S.,  "Measurements  in  Turbulent  Boundary  Layers,"  M.S.  Thesis, 
Dept.  Aerospace  and  Ocean  Engrg.,  VPI&SU,  in  progress,  1985. 

23.  McGrath,  B.,  "Some  Features  of  Surface  Pressure  Fluctuations  in  Zero 
and  Favorable  Pressure  Gradients,"  M.S.  thesis.  Dept.  Aerospace  and 
Ocean  Engrg.,  VPI&SU,  in  progress,  1985. 

24.  Chew,  Y.-T.  and  Simpson,  R.L.,  "A  Simple  Data  Reduction  Method  for 
Triple  Sensors  Hot-Wire  Anemometers  in  Three-Dimensional  Flow," 
manuscript  in  preparation,  1985. 


PUBLICATIONS 


Pelletier,  D.  and  Schetz,  J.A.,  "A  Navier-Stokes  Calculation  of  3D 
Turbulent  Flow  Near  a  Propeller  in  a  Shear  Flow,"  AIAA  Paper  85-0365, 
January  1985. 

Pelletier,  D.  and  Schetz,  J.A.,  "A  Turbulence  Model  for  Finite  Element 
Simulation  of  3-D  Turbulent  Flows  Near  Propellers  and  Windmills," 
Fourth  Int .  Conf.  Numerical  Methods  in  Laminar  and  Turbulent  Flow, 
Swansea,  Wales,  July  1985. 


p 

R 

5 

5 

U 

R 

E 

C 

Q 

E 

F 


0.3 

1 


LEOE'iO:  Z  -  X/0=0.025 

.  X/D=0. 100 

-  X/D=0.23 

-  X/0=a. SIS 

- X/D  =  0.  9311 

-  X/0=1 . 5u 

-  X/n  =  2.'4? 

Figure  15:  Shear  Flow:  Pressure 
Profiles 


CONCLUSIONS 


This  work  has  developed  a  general- 
purpose,  advanced  computational  techni¬ 
que  for  thrse-di.mensional  turbulent  free 
shear  flows  past  prcpellsrs.  The  main 
results  can  be  summarized  as: 

1.  predicticns  for  uniform  and  shear 
flows  past  a  propeller  show  ex¬ 
cellent  agreement  with  the  e.xper- 
iments  for  the  axial  and  swirl 
component  of  the  velocity, 

2.  the  simple  ourbulence  model  pro¬ 
duced  vary  accurate  swirl  predic¬ 
tions  for  both  the  uniform  and 
shear  flow  cases.  All  previous 
numerical  work  on  such  flows  pro¬ 
vided  poor  prediccions  of  the 
swirl, 

3.  the  pressure  predictions  present¬ 
ed  here  for  viscous  turbulent 
flows  through  a  propeller  are  the 
first  known  to  the  authors. 

4.  the  numerical  simulations  show 
that  the  propeller  e.xerts  a 
strong  upstream  influence  on  the 
axial  and  radial  velocity  fields. 


and  on  the  pressure  field.  There 
is  negligible  swirl  upstream  of 
the  propeller, 

5.  the  resulting  _aierkin  finite 
element  algorithm  is  robust  and 
stable,  mostly  due  to  the  full 
coupling  of  the  three  momentum 
equations  and  the  continuity 
equation.  It  captured  the  com¬ 
plex  features  of  the  shear  flow. 
No  gradual  i.ntroduction  of  the 
thrust  and  torque  as  the  itera¬ 
tion  proceeds  is  required  as  op¬ 
posed  to  previous  work. 
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peller  on  the  pressure  field  is 
illustrated  in  Figure  (10).  The  sudden, 
almost  discontinuous,  pressure  increase, 
characteristic  of  actuator  disks,  is  ev¬ 
ident.  Figure  (11)  shows  that  there  is 
negligible  swirl  upstream  of  the  pro¬ 
peller  as  predicted  by  vortex  theory. 
The  jump  in  the  swirl  across  the  pro¬ 
peller  is  a  behavior  akin  to  the  shock 
wave  in  compressible  flows.  The  small 
wiggles  upstre2un  of  the  propeller  are 
typical  of  the' use  of  Galerkin  and  ce.n- 
tral  difference  sche.mes  for  problems 
with  sharp  fronts.  The  near  perfect  ax¬ 
ial  symmetry  of  the  swirl  distribution 
is  observed.  The  slipstream  contraction 
is  clearly  seen  in  Figure  (12) 

Figure  (13)  presents  axial  velocity 
profiles  along  the  y  a.xis  at  various 
stations  downstream  of  the  propeller. 
The  streamwise  acceleration  of  the  fluid 
downstream  of  the  propeller  due  to  pres¬ 
sure  relaxation  can  be  seen.  The  kinks 
near  the  tip  of  the  propeller  in  the 
first  two  profiles  of  Figure  (13)  are 
characteristic  of  flows  past  propellers 
?,ef.[151.  The  asymmetries  due  to  the 
shear  in  the  approach  flow  can  be  seen 
most  clearly  in  the  shape  of  the  peaks 
of  axial  velocity.  The  swirl  profiles 


at  various  axial  stations  are  shown  in 
Figure  (14).  The  swirl  level  decays  very 
slowly.  The  radial  turbulent  diffusion 
of  momentum,  characterized  by  flatter 
and  wider  profiles,  is  clearly  seen  . 
Pressure  profiles  are  plotted  in  Figure 
(15)  illustrating  the  decay  downstream 
of  the  propeller.  The  negative  dimples 
in  the  last  three  curves  are  due  to  the 
zero  axial  traction  applied  at  the  out¬ 
flow.  The  use  of  this  natural  boundary 
condition  sets  the  pressure  level  at  the 
outflow  boundary.  The  axial  traction  is 
given  by: 

t  =  -  p  +  2u_U,  =  0 
X  T  X 

where  t^  is  the  axial  traction  and  P  the 
pressure.  The  zero  traction  defines  the 
pressure  level  as: 

P  =  2u^U,^ 

Since  the  axial  velocity  is  still  decay¬ 
ing,  U,^  is  negative  and  the  pressure  at 
the  exit  takes  a  small  negative  value. 
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figure  14:  Shear  Flow:  Swirl  Profiles 


Figure  13;  Shear  Flow.  Axial  Velocity 
Profiles 
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Figure  5:  Uniform  Flow:  Pressure 


K,  =  T/fpQN^D*'^)  =  0.153 

Kg  =  Q/(Pq.N^D*^)  =  0.0266 

.^t  the  inflow  and  on  the  freestream 
shell  boundary  ,  U  is  set  equal  to  the 
approach  shear  flow.  All  other  boundary 
conditions  are  identical  to  those  used 
for  the  uniform  flow  simulation. 

To  obtain  an  initial  guess  of  the 
velocity  field,  a  2-D,  axisymnetric  so¬ 
lution  without  swirl  was  produced  for  a 
uniform  flow  and  the  velocity  excess  was 
extracted,  rotated  on  the  3-D  grid,  and 
added  to  the  approach  shear  flow.  The 
system  of  nonlinear  equations  was  then 
solved  in  2  quasi-Mewton  iterations  re¬ 
quiring  a  total  of  150  .minutes  of  CPU  on 
I3M-30ai.  No  gradual  introduction  of 
thrust  and  torque  in  the  course  of  the 
iteration  was  required  as  opposed  to  the 
work  of  Refs. [4-6] . 

Figure  (6)  compares  the  experimental 
and  predicted  axial  velocity  profiles 
along  the  y  axis.  The  agreement  is  ex¬ 
cellent  at  x/D  =  0.025  and  good  at  x/D  = 
0.23.  On  the  low  speed  side  of  the  pro¬ 
peller,  less  fluid  flows  through  the 
propeller,  and  hence,  it  is  subject  to 
more  body  force  per  unit  mass  flow.  This 
results  in  higher  acceleration  of  the 
fluid  and  a  larger  relative  speed 
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Figure  6:  Shear  Flow:  .Axial  Velocity 


increase  on  the  positive  side  of  the  y 

axis  compared  to  the  negative  side. 

Figure  (7)  presents  the  predicted  and 

experimental  swirl  profiles  along  the  z 

axis.  Agreement  between  predictions  a.nd 

experi.ments  is  excellent.  Figure  (3) 

gives  the  predicted  and  e.xperimentai 

pressure  profiles  along  the  y  axis  in 

the  form  of  a  pressure  coefficient  C  = 

? 

2 

2(P-Pg)/(pUg)  .  The  experimental  pro¬ 
files  have  shapes  similar  to  those  of 
the  predictions,  but  the  amplitudes  do 
not  quite  agree  at  x/D=0.025. 


Ficfures  (9)  to  (12)  present  predicted 
axial  distributions  of  the  axial  compo¬ 
nent  of  the  velocity,  pressure,  swirl, 
and  radial  velocity.  These  distribu¬ 
tions  were  taken  at  a  radius  equal  to 
0.35R  where  the  maxi.mum  thrust  occurs. 
Figure  (9)  clearly  shows  the  strong 
streamwise  acceleration  of  the  fluid  up¬ 
stream  and  downstrea.m  for  the  propeller. 
The  different  variations  in  U  at  four 
peripheral  location  are  due  to  the  shear 
flow.  The  two  curves  in  the  vertical 
plane  display  non-symmetrical  behavior. 
This  is  caused  by  the  swirl  that  moves 
high  velocity  fluid  into  regions  of  low¬ 
er  velocity  and  vice  versa .  The  signi¬ 
ficant  upstream  influence  of  the  pro- 
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A  note  of  caution  is  required  before 
any  comparison  with  experiments  is  at¬ 
tempted.  While  the  finite  element  simu¬ 
lations  can  accurately  represent  a  free 
running  propeller,  the  experiments  must 
use  a  shaft  to  support  and  drive  the 
propeller.  A  body  must  be  placed  close 
to  and  downstream  of  the  propeller  to 
house  the  drive  train,  see  Ref. [16). 
Xhe  housing  begins  at  0.23  propeller 
diameters  downstream  of  the  disk.  The 
presence  of  the  housing  will  very  likely 
affect  the  experi-mental  results  from  the 
station  located  at  x/D  =  0.23,  since 
this  station  sits  right  on  the  body. 

Figure  (3)  compares  the  experimental 
and  predicted  axial  velocity  profiles 
along  the  y  axis.  The  agreement  is  ex¬ 
cellent  at  x/0  =  0.025  and  good  at  x/D  = 
0.23. 

Figure  (4)  presents  the  predicted  and 
experimental  swirl  profiles  along  the  y 
axis,  rt  can  be  seen  that  the  assumed 
radial  distribution  of  the  swirling  body 
force  is  reasonable,  but  it  has  slightly 
too  sharp  a  pea.k  located  too  close  to 
the  propeller  tip.  Agreement  between 
predictions  and  experiments  is  excellent 
and  constitutes  a  major  improvement  over 
those  of  Refs. [4-61. 

Figure  (5)  presents  the  predicted  and 

experi.mental  pressure  profiles  along  the 

y  axis.  The  pressure  is  given  in  the 

form  of  a  pressure  coefficient,  C  = 

P 

2(?-?q)/(pUq) .  These  are  the  first 
viscous  pressure  predictions  for  this 
type  of  flow  known  to  the  authors.  All 
analysis  methods  based  on  the  boundary 
layer  equations  require  the  specifica¬ 
tion  of  the  pressure  as  input.  The 
agreement  with  the  experiments  at 
x/D=0.02S  is  good.  The  pressure  pred¬ 
icted  at  x/D=0.23  is  in  qualitative 
agreement  with  the  theory  of  free  run¬ 
ning  propellers.  The  experi.mental  pro¬ 
files  at  x/0=0.23  have  negative  values 
typical  of  the  flow  over  an  obstacle 
such  as  the  drive  housing. 
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Figure  3:  Uniform  Flow:  Axial  Velocity 
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SHEAR  FLOW  PAST  A  PROPELLER 

This  si.mulation  of  a  free  running 
propeller  in  a  shear  flow  corresponds  to 
the  experiments  of  Kotj  [161.  The  char¬ 
acteristics  of  the  flow  are  as  follows; 

U*j^j  =  8.52  -  0.67y  m/s 
=  1  -  0.1588  y 

p*  =  1.177  kg/m? 

D*  =  0.492  m 

N  =  1150  rpm. 
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four  curves  result  in  clustered  distri¬ 
butions  well  within  the  margin  of  uncer¬ 
tainty  for  this  problem.  An  initial  va- 
lue  of  5*10  ■  was  judged  adequate  and 
used  for  all  subsequent  simulations. 
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Figure  1:  Effect  of  Initial  Value  of 
Eddy  Viscosity 

CCHPITTATIONAL  SAVINGS  ARISING  FROM  THE 
USE  OF  THE  PENALTY  METHOD . 

To  illustrate  the  computational  sav¬ 
ings  afforded  by  the  use  of  the  penalty 
method,  we  present  the  characteristics 
of  the  global  stiffness  matrix  of  the 
3-D  propeller  simulations  for  the  penal¬ 
ty  and  velocity-pressure  (mixed)  formu¬ 
lations.  The  mesh  has  a  total  of  3768 
nodes  and  3588  8-noded  bric)«  elements. 


Penalty 

Mixed 

NEQ 

10  557 

14  145 

NEM 

10  997  103 

19  434  905 

MINE 

521 

687 

MAXB 

886 

- 

1  186 

where  NEQ  is  the  nu.mber  of  equations, 
NE.M  the  number  of  matrix  coefficients 
stored,  MIN3  the  average  bandwidth  of 
the  matrix,  and  MAXB  the  maximum  band¬ 


width.  The  mixed  formulation  has  34% 
more  unknowns,  a  bandwidth  32%  greater, 
and  stores  77%  more  matrix  coefficients. 
The  Gaussian  L-U  factorization  time  is  2 
and  5  hours  CPU  for  the  penalty  and  mix¬ 
ed  methods  repectively  (on  IBM-3081,  us¬ 
ing  FORTRAN-H  extended  with  optimization 
level-1).  Clearly,  the  penalty  method 
offers  substantial  savings. 

UNIFORM  FLOW  PAST  A  PROPELLER 

The  characteristics  of  the  flow  are 
as  follows; 

U*  =  8.52  m/s 

p*  =  1.177  )tg/m^ 

D*  =  0.492  m 

N  =  1150  rpm 

K.J,  =  T/(p*N-^D*^)  =  0.150 

Kg  =  Q/(pqN^D*^)  =  0.0279 

The  grid  is  a  cylinder  with  24  nodes  in 
the  axial  direction,  14  in  the  radial 
direction  ,  and  12  in  the  peripheral  di¬ 
rection.  The  inflow  boundary  is  located 
2  diameters  upstream  of  the  propeller, 
see  Fig. (2).  The  outflow  plane  stands 
at  3  diameters  downstream  of  the  dis)t, 
while  the  freestream  boundary  is  a  cyl¬ 
indrical  shell  of  radius  1.2  diameters. 
At  the  inflow,  U  is  set  to  unit  value 
and  V  and  W  are  set  to  zero.  On  the 
freestream  boundary,  U  is  set  to  1  and 
the  y  and  z  tractions  to  zero.  The  chree 
components  of  the  traction  vector  are 
set  to  zero  at  the  outflow  boundary. 


Figure  2:  Grid  for  Propeller  Simu'.ation 


RESULTS  AND  DISCUSSION 

The  turbulence  model  was  programmed 
and  interfaced  to  the  general  purpose 
finite  element  fluid  dynamics  program 
elDAP  (15).  The  resulting  code  can  simu¬ 
late  planar,  axisymmetric,  and  three-di¬ 
mensional  turbulent  free  shear  flows. 


of  the  propeller  and  iU  is  the  charac¬ 
teristic  velocity  excess  determined  from 
a  simplified  inviscid,  one-dimensional 
analysis  of  the  propeller  (131: 

AU  =  (  2T  /(p  A)  -  \  ^  -  Uq 

The  following  data  applies  to  the  uni¬ 
form  flow  case  of  Ref. [161: 

R*  =  0.246  m 
U*  =3.52  m/s 


NOMDI MENS I ONAL I 2AT ION 


=  2.914  N 
=  0.19  m^ 

=  1.177  kg/m^ 


All  calculations  were  performed  with 
a  nondimensicnal  form  of  the  equations 
of  motion  for  the  conditions  of 


conditions 


Ref. [161.  Reference  values  are  selected 

*  t 

for  the  velocity  (U^),  length  (Lq), 

* 

and  pressure  (Pq)  to  obtain  the  fol¬ 
lowing  dimensionless  variables  (the  star 

denotes  a  dimensional  variable): 

*  * 

Xi  =  x./Lq 

U,  =  u*/u;^ 

?  =  (?*-?*)/(puJ) 

U«  =  u.j*/(pU*L*)  = 

f  =  fV(uJ/L;) 

p  =  p/pq  =  : 

INITIAL  VALUE  FOR  THE  EDDY  VISCOSITY 


The  absence  of  solid  walls  or  aftbody 
upstream  of  the  propeller  in  the  pre¬ 
sent,  idealised  flow  problem  results  in 
an  undisturbed  "inviscid",  uniform  or 
shear  approach  flow  making  difficult  the 
determination  of  the  initial  condition 
for  the  turbulence  model.  The  flow 
through  and  past  a  propeller  has  many  of 
the  characteristics  of  a  turbulent  jet 
[17].  Hence,  an  indication  of  the  mag¬ 
nitude  of  may  be  obtained  from  the 

A 

standard  jet  formula  for  the  eddy  vis¬ 
cosity  [171: 


V.J  =  0.025  r^^2 


where  ‘-^^2  halfwidth  of  the  layer 

and  can  be  taken  as  radius 


AU  =1.41  m/s 

A  rough  estimate  of  the  eddy  viscosity 
downstream  of  the  propeller  is  in  dimen¬ 
sionless  form: 


=  Uj/(DUq)  =  0.002 


There  remains  to  select  an  initial 
value  for  the  eddy  viscosity.  Assuming 
that  the  Prandtl-Kol.mogorov  relationship 
holds  we  have: 

v_  =0.2  kT'  R 
^init  “ 

where  c=0.2  (the  value  for  jets)  and  k^ 
is  the  turbulence  intensity  upstream  of 
the  propeller.  The  corresponding  di.men- 
sionless  form  is 


T.  . 
xnit 


P*  /(UqD)  =  0.1  kJ^VUc 


From  the  experiments  of  Ref. [15]  we  have 

and  an  esti.mate  of  the  initial  eddy  vis¬ 
cosity  is  finally 

V,.  =  C.COl 


To  assess  the  validity  of  this  analy¬ 
sis  and  to  investigate  the  effect  of  the 
choice  of  the  initial  value  on  the  eddy 
viscosity  distribution,  simulations  were 
performed  for  the  uniform  flow  of 
Ref. [161.  The  calculations  were  done 
for  an  axisymmetric  flow  without  swirl. 

Figure  (1)  shows  that  low  values  of 
in  the  range  lO"^  to  lO'^,  produce 


in  the  range  lO"^  to  lO'^,  produce 
similar  results  with  a  peak  eddy  viscos¬ 
ity  of  0.0025,  a  value  in  good  agreement 
with  our  rough  estimate  of  0.002.  High¬ 
er  initial  values  produce  distributions 
that  are  probably  too  high.  The  bottom 
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where  are  the  residuals,  a  measure  of 
the  quality  of  the  approximate  solution 
used.  The  Galerkin  method  reduces  this 
error  to  zero,  in  a  weighted  sense,  by 
making  the  residuals  orthogonal  to  some 
set  of  functions.  The  weighting  func¬ 
tions,  must  satisfy  the  continuity 
equation.  The  continuity  equation 
constraint  is  enforced  by  use  of  a  pe¬ 
nalty  method.  The  resulting  weak  formu¬ 
lation  is  given  by  [7,11]: 


For  reasons  of  computational  economy, 
we  use  the  trilinear  velocity,  constant 
pressure  element  [11].  This  element, 
under  certain  circumstances  may  suffer 
from  spurious  pressure  solutions.  Our 
experience  with  this  element  indicates 
that  it  is  a  fairly  reliable  element  for 
the  problems  treated  here.  No  spurious 
pressure  solutions  were  observed. 

GENERAL  SOLUTION  PROCEDURE 


1 5V^ ,IJ.  U,  ^  .  *5U .  _  ( U .  _  ^  *V.  _  .  )  +«U .  pf .  ! dE 

^D  '^i,i®^i,i  =  ■fs  *“i^i 


Here  X  is  a  positive  nu-mber  whose  value 
affects  the  accuracy  of  the  solution, 
while  5U^  may  be  interpreted  as  a  virtu¬ 
al  velocity,  and  aq.(lO)  is  similar  to 
the  principle  of  virtual  work  of  solid 
The  divergence  theorem  was  applied  to 
the  pressure  and  viscous  terms  to  reduce 
the  differentiability  requirement  on  the 
velocity  by  transfering  some  of  the  der¬ 
ivative  from  U^  to  5U^,  and  to  introduce 
the  natural  boundary  conditions  involv¬ 


ing  the  surface  tractions  or  forces 

ti  =  (-?  * 


The  pressure  does  not  appear  explicitly 
in  eq.(lO);  it  appears  only  implicitly 
through  the  surface  integrals  on  the 
right  hand  side  of  eq.(lO).  The  only 
unknowns  are  the  velocity  components. 
This  results  in  significant  computation¬ 
al  savings.  Once  the  velocity  field  is 
obtained  by  solving  eg. (10),  the  pres¬ 
sure  may  be  computed  in  a  post  process¬ 
ing  step  from  the  following  relationship 
[7,111: 


P  =  -XU^  ^  (11) 

Convergence  of  the  solution  of  the 
penalty  formulation  to  the  true  solu¬ 
tion,  as  X-»,  can  be  proved  [12].  A  va- 

£  O 

lue  of  X  of  10  to  10  usually  proves 
effective . 


The  formulation  results  in  a  system 
of  nonlinear  algebraic  equations  of  the 
form: 

(K(U,u^(U))]  U  =  F 

where  K  is  the  global  system  matrix,  U 
is  the  global  vector  of  unknowns,  and  F 
is  the  global  vector  of  body  forces  and 
boundary  conditions. 

This  system  of  equations  may  be 
solved  by  a  combination  of  methods: 
successive  substitution,  Newton-Raphson, 
and  quasi-Newton  methods  [13].  The  qua¬ 
si-Newton  procedure  converges  as  fast  as 
the  full  Newton  method,  but  the  cost  of 
one  quasi-Newton  iteration  is  aibout  10 
to  15%  cf  that  of  a  Newton  iteration. 
The  following  iteration  strategy  was 
found  suitable  throughout  this  study: 
starting  from  a  first  guess  for  the  vel¬ 
ocity  field,  the  quasi-Newton  method  is 
used  to  Iterate  until  convergence.  The 
eddy  viscosity  distribution  is  updated 
at  each  iteration  by  solving  eq.(6)  us¬ 
ing  the  current  %’elocity  distribution  to 
evaluate  the  coefficients  of  the  ordi¬ 
nary  differential  equation.  The  linear¬ 
ized  system  of  equations  is  very  large 
and  is  solved  by  direct  Gaussian  decom¬ 
position  in  a  compacted,  skyline,  out  of 
core  solver  [ 14 j . 
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dv^/dx  =  O.SCY^-v^dY^/dx-v^*  <S) 

where 

Y^  =  11^  UlV(C^L)  dA 

Y3  =  //^  ra2/(cjL)  dA 

It  now  remains  to  specify  the  length 
scale  L.  Following  Ref. [91  we  adopt: 

^  =  ^<^1/2  2i/25/(^!/2  *  2^/2)'"''* 

with  a  =  1.89.  This  ensure  correct 

dependence  for  the  limiting  cases  of 
planar  and  axisymmetric  flows.  Equatio 
(6)  is  an  initial  value  problem  for  the 
kinematic  eddy  viscosity  The  only 

data  required  is  an  initial  condition 
for  at  the  upstream  inflow  boundary, 
a  value  that  can  be  easily  estimated  for 
most  flows.  One  significant  advantage  of 
the  present  .model  is  the  absence  of  the 
diffusion  of  TKE,  a  most  difficult  term 
to  model. 

The  present  model  was  calibrated  on 
the  simple,  well  documented  turbulent 
flow  problem  of  the  far  field  of  a  round 
]et  issuing  into  still  surroundings. 
This  flow  has  an  analytical  solution 
;<SR11.)  for  the  velocity  field,  and  the 
value  of  the  eddy  viscosity  is  known  and 
constant.  The  integrals  are  evaluated 
exactly  and  with  the  assu.mption  of  tur¬ 
bulence  equilibrium  the  constants  are 
determined  to  be  [ 7  ]  ; 

83  =  C.519  Cj  =  0.154 

MODELING  THE  PROPELLER 

The  propeller  is  modeled  by  a  disk  of 
radius  equal  to  the  propeller  radius  and 
of  thickness  Ax,  roughly  equal  to  the 
physical  thickness  of  the  propeller.  The 
thrust  and  torque  are  allowed  to  vary 
radially  but  are  constant  in  the  tangen¬ 
tial  direction.  Little  is  known  about 
the  detailed  radial  thrust  and  torque 
distributions  of  a  given  propeller. 
Generally,  one  only  knows  the  global  va¬ 
lues  of  total  thrust  and  torque.  For 


simplicity  we  use  a  trapezoidal  distri¬ 
bution  given  by: 

t(r)  =0  r  in  '0,rll 

t(r)  =  t^(r-r3)/(r2-rj^)  r  in 

t(r)  =  tjjj  r  in  [r2,r3]  (7) 

t(r)  =  tj^(R-r)/(R-r3)  r  in  [r3,Rl 

t(r)  =  0  r>R 

where  t  is  the  maximum  value  of  the 
n 

thrust  and  R  is  the  radius  of  the  pro¬ 
peller.  Values  of  rj^,  r2/  and  r3  were 
set  to  0.25R,  0.7R,  and  0.85R  respec¬ 

tively.  These  values  result  in  distri¬ 
butions  similar  to  those  of  Refs. [4-51. 
The  same  form  is  adopted  for  the  dictri- 
bution  of  the  force  s  producing  the 
swirl.  Its  maximum  value  is  denoted  by 
s^.  The  upstream  face  of  the  propeller 
disk  is  located  at  x  =  -  ix,  while  the 
backface  of  the  propeller  is  located  at 
X  =  0.0.  These  distributions  are  inte¬ 
grated  to  yield  the  global  thrust  and 
torque  of  the  propeller  : 


T  =  0.3075 

*  2 

II  * 

ft 

B 

(S) 

Q  =  0.2218 

*  2 

n  * 

R^  *  s 

m 

(9) 

Given  values 

of 

the 

global  thrust 

and 

torque,  the  values  of  t  and  3_  are  det- 

mm 

ermined  form  eq3.(8)  and  (9). 


SOLUTION  ALGORITHM 


PENALTY  FORMULATION 

Details  of  the  weak  Galerkin  formula¬ 
tion  are  readily  available  for  penalty 
function  formulations  [11).  Thus,  only 
an  outline  of  the  technique  is  presented 
here . 

Substitution  of  an  appro.ximate  solu¬ 
tion  (U  ,P  )  into  eqs.(l)  and  (2) 
yields  a  set  of  residual  equations  of 
the  form: 

*  * 

Momentum  :  ,?  )  =  R. 

Continuity  :  f2(U*)  =  R2 
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MODELING  OF  THE  PROBLEM 


^i,i  =  ° 

^  <2) 

where  is  the  turbulent  eddy  viscosi¬ 
ty,  and  are  body  forces  representing 
the  effects  of  the  propeller  disk. 

TURBULENCE  MODEL 

Following  Refs. (4-61  a  simple  turbu¬ 
lence  model  is  used,  since  little  e.'cper- 
imental  data  is  available  for  this  type 
of  problem  on  which  to  base  more  elabo¬ 
rate  models.  An  integrated  form  of  the 
TKE  equation  is  used  (details  of  the 
algebraic  manipulations  are  contained  in 
Ref . [71 ) : 


where  L  is  a  characteristic  length  scale 
of  the  shear  layer,  and  ^2  and  C2  are 
constants.  Physically,  all  integrands  in 
aq.(3)  are  sufficiently  well  behaved  to 
ensure  boundedness  of  the  integrals.  The 
assumption  that  the  eddy  viscosity  is 
constant  over  the  cross  section  of  the 
flow  is  a  well  documented  behavior  for 
free  shear  flows.  In  fact,  the  eddy  vis¬ 
cosity  is  constant  over  most  of  the 
shear  layer  and  decays  to  zero  only  as 
the  radial  distance,  r,  from  the  x-axis 
goes  to  infinity.  In  order  to  preserve 
this  behavior  and  ensure  boundedness  of 
the  integrals  in  the  modelled  integrated 
equation,  we  introduce  a  distribution 
function,  T,  representing  the  distribu¬ 
tion  of  Wj  across  the  layer.  It  is 
closely  related  to  the  actual  intermit- 
tency  of  the  turbulent  flow.  The  fol¬ 
lowing  form,  obtained  from  a  nonlinear 
least  squares  fie  to  intermittency  data 
for  a  turbulent  boundary  layer  over  a 
flat  plate  and  turbulent  plane  and  round 
jets  was  found  satisfactory: 

where 

Xy  =  0.5  (  1  -  arf(1.98Y/Y^^2'2-'‘2)  > 

*  0.5  (  1  -  erf(1.98Z/f*^2"3-‘*2)  ) 


d(/;^pUk  dA)/dx  = 

-  dA  (3) 

In  order  to  close  the  .model,  the  TKE 
must  be  related  to  the  eddy  viscosity 
u._=pv.^.  This  is  done  through  the  use  of 
the  Prandtl-Kolmogorov  relationship  (8] 

=  C2  k^/^  L  (4) 

Assuming  that  laminar  dissipation  is  im¬ 
portant  only  in  the  inertial  sub-range, 
dimensional  consistency  dictates  that 
the  viscous  dissipation  be  modeled  as 
(81  : 

f  =  a2  k^/^  /  L  (5) 


and  ^1/2  half-widths  of 
the  shear  layer  along  the  y  and  z  axes 
respectively.  The  half-widths  are  de¬ 
fined  by  the  points  on  the  axis  where  : 


In 


(U  -  ■ 

this  expression, 


u  ' 

inf 'max 

'^'inf 


0.5 

the  free 


stream  approach  flow  velocity  which  can 


be  non-uniform. 


It  should  be  noted  that  the  particu¬ 
lar  form  of  the  distribution  function 
used  is  not  critical.  Choosing  diffe¬ 
rent  algebraic  forms  will  simply  result 
in  slightly  different  values  of  the 
constant  a2  and  C2. 

Upon  substitution  of  eq.(4)  and 
eq.(S)  into  eq.(3),  the  modeled  inte¬ 
grated  TKE  equation  takes  the  form: 
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INTRODUCTION 


Detailed  analysis  of  the  3-D,  turbu¬ 
lent  flow  produced  by  a  propeller/body 
combination  is  of  interest  in  a  number 
of  practical  applications.  Examples  in¬ 
clude  the  influence  of  the  propeller  on 
the  body  pressure  distribution,  the 
prediction  of  the  near-wake  profiles, 
cyclic  loading  to  produce  vibrations  and 
the  influence  of  a  downstream  surface  *on 
a  propeller  and  vice  versa,  to  name  a 
few.  Until  recently,  only  approximate 
treatments  were  available  involving  one 
or  more  of  the  following  restrictive  as¬ 
sumptions.-  the  flow  was  assumed  invis- 
cid;  the  propeller  was  represented  as  an 
actuator  disk  with  constant  thrust;  the 
flow  was  taken  as  laminar  and/or  the  ef¬ 
fects  of  the  propeller  on  the  flowfield 
were  assumed  small  enough  to  permit  li¬ 
nearization  of  the  equations  of  motion. 
Refs. (1-3)  are  representative. 

The  numerical  solution  procedure 
reported  in  Refs.  (4-6)  had  as  its  goal 
the  development  of  a  realistic  treat¬ 
ment,  holding  si-mplifying  assumptions 
and  approxi.mations  to  a  minimum.  The 
work  was  based  on  the  full,  axisymme- 
tric,  mean  (in  the  turbulence  sense), 
unsteady  Navier-Stokes  equations.  In 
order  to  place  some  bounds  on  the  scope 
of  the  effort  at  that  ti.me,  however, 
some  simplifications  were  necessary. 
The  first  was  the  assumption  of  an  actu¬ 
ator  disk  model  for  the  propeller,  alt¬ 
hough  arbitrary  radial  variations  of 
thrust  were  allowed.  Second,  the  flow 
was  taker,  as  axisymmetric .  Third,  turbu¬ 
lent  transport  processes  were  described 
by  an  integrated,  turbulence-kinetic-en¬ 
ergy  (TKE)  model,  which  was  used  to 
predict  an  eddy  viscosity  distribution. 
The  unsteady  equations  of  motions  were 
cast  in  ter.ms  of  a  stream  function,  one 
vorticity  component,  and  the  peripheral 
velocity.  Comparison  of  the  predictions 
of  this  procedure  with  laboratory  data 
for  an  axisymmetric  flow,  showed  good 


agreement  for  the  axial  velocity  compo¬ 
nent.  The  swirl  velocity  component  pred¬ 
ictions  were  consistently  too  low. 

From  the  point  of  view  of  direct  ap¬ 
plications  to  problems  of  practical  in¬ 
terest,  the  biggest  limitation  to  the 
analysis  described  above  is  its  restric¬ 
tion  to  two-dimensional,  axisymmetric 
flows.  Actual  propeller-driven  vehicles 
have  either  a  three-dimensional  body 
near  the  propeller  and/or  appendages 
that  render  the  flow  three-di.mensional . 
For  many  cases,  the  assumption  of  an  ac¬ 
tuator  disk  representation  of  the  pro¬ 
peller  remains  appropriate  while  the 
restriction  to  axisymmetric  flow  does 
not.  Furthermore,  the  jump  to  trying  to 
treat  the  three-dimensional,  cyclicaly 
unsteady  problem  with  individual  blades 
and  their  thin  boundary  layers  and  wakes 
is  too  great  to  be  attempted  in  one  step 
at  this  time. 

The  present  work  is  part  of  a  step- 
by-step  approach  to  the  development  of  a 
computational  method  for  analysis  of 
propeller  flowfields  with  three-dimen¬ 
sional  inflows  and  the  affects  of  ap¬ 
pendages  in  situations  of  increasing 
complexity  and  practical  realism.  The 
fully  elliptic,  three-dimensional,  ti.me 
averaged,  steady  state,  primitive  varia¬ 
bles  Mavier-Stokes  equations  are  solved 
by  a  Penalty  Finite  Element  ■•lethod  be¬ 
cause  of  its  ease  of  handling 
complicated  geo.metries  and  its  simplicty 
for  implementing  a  variety  of  boundary 
conditions.  Turbulence  modeling  is  done 
through  a  generalization  of  the  inte¬ 
grated  TKE  model  of  Refs. (4-5]. 
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NOMENCLATURE 


Abstract 


A  numerical  procedure  based  on  the 
primitive  variable  Navier-Stokes 
equations  is  applied  to  the  simulation 
of  the  three-dimensional  flow  near  a 
propeller  in  a  shear  flow.  The 
Navier-Stokes  equations  are  solved  by  a 
Penalty  Finite  Element  Method.  The 
propeller  is  .modelled  as  an  actuator 
disk,  and  the  direct  simulation  of  a 
given  propeller  is  considered  in  detail. 
Turbulent  transport  is  modeled  by  an 
integrated  Turbulent  Kinetic  Energy 
equation.  This  approach  results  in  a 
robust  numerical  algorithm.  Detailed 
comparison  with  wind  tunnel  measurements 

show  good  prediction  of  velocity  and 
pressure.  The  high  accuracy  of  the  swirl 
prediction  is  a  major  i.T.provement  over 
previous  analyses. 


A  . . .  cross  section  of  the  flowfield  and 
area  of  the  propeller 
•  •  constant  in  the  turbulence  .model 
C2  . .  constant  in  the  turbulence  model 
D  . . .  domain  of  solution  and 
propeller  diameter 
K  . . .  turbulence  kinetic  energy 
Kq  . .  torque  coefficient 
. .  thrust  coefficient 
N  .  ...  rotational  speed  of  the  propeller 
n^  . .  unit  outward  normal 
P  ...  pressure 

Q  ...  torque  of  the  propeller 
^2,/2  half -width  of  shear  layer 
Re  . .  Reynolds  number 

s  ...  radial  distribution  of  swirl  force 
S  ...  boundary  of  the  domain  D 
t  ...  distribution  of  thrust 
t^  . .  surface  traction  vector 
T  ...  thrust  of  the  propeller 

U^  . .  velocity  vector 

U^^j  upstream  approach  flow 

p  ...  density  of  the  fluid 

t  ...  turbulence  viscous  dissipation 

i  ...  penalty  parameter 
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